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We present a simple and efficient method to compute the superhorizon evolution of the spectral 
index in multi-field inflationary models, using transport equation techniques. We illustrate the 
evolution of ris with time for various interesting potentials. 
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I. INTRODUCTION 

The imminent arrival of precise microwave background 
data from Planck makes this a promising time to put 
constraints on cosmological models, and in particular the 
inflationary paradigm. To do so, we need accurate cal- 
culations of the observables produced by an inflationary 
phase, of which we expect the spectral index of the den- 
sity perturbation, ris, to be among the most precisely 
measured. Planck may determine rig with an error of 
perhaps a few times 10"'^ [T]. 

Inflationary models motivated by high-energy physics 
often invoke many light scalar fields. The resulting multi- 
field dynamics can cause the curvature perturbation syn- 
thesized at horizon crossing to evolve, as power is trans- 
ferred from decaying isocurvature modes. During this 
time the spectral index will typically vary. Accurately 
tracking its trajectory is an important challenge, and var- 
ious techniques to carry out the computation have been 
proposed. 

Statistical properties of the inflationary density per- 
turbations can be calculated using the separate universe 
approximation. A simple way to implement this approx- 
imation is to Taylor expand in powers of the initial con- 
ditions When applied to the curvature perturbation 
this is sometimes called the "(5iV formalism" . Long ago, 
Sasaki & Stewart used this method to provide a very 
general formula for the spectral index valid for an 
arbitrary choice of field space. This formula has proved 
entirely satisfactory for analytic calculations. But our 
ambition to constrain increasingly complex models means 
that numerical work is often needed, particularly for the 
sophisticated examples motivated by high-energy the- 
ory in E] . For these CcLS6S, 3j direct implementation of the 
5N formula is less attractive, because it is not framed in 
terms of ordinary differential equations but rather "vari- 
ational" derivatives which require numerical integration 
of high accuracy. For a discussion, see Ref. [6] . 

An alternative approach was developed by Gordon et 
al. and by Nibbelink & van Tent [7]. Their formalism 
is equivalent to the Taylor expansion method, although 
more complicated because it requires decomposition in 
terms of a Frenet basis for the inflationary trajectory 
in field space. Recently, Peterson & Tegmark used this 
method to compute the spectral index JF. But it would 
be desirable to have a simple system of ordinary differ- 
ential equations, expressed directly in terms of the po- 
tential, which do not invoke these complications. In this 



short paper we provide such a system. 

To achieve this, we explore the evolution of correla- 
tion functions in a multi-field framework using transport 
techniques [SI IH] • We obtain a transport equation which 
propagates the spectral index from its value at horizon 
crossing up to the end of inflation. The result is tech- 
nically equivalent to those of Stewart & Sasaki and Pe- 
terson & Tegmark, but we believe it has advantages of 
simplicity and numerical implementation. In addition, 
it is easy to study the entire history of the spectral in- 
dex during inflation, which can be correlated with other 
observables such as the amplitude of the local-mode bis- 



pectrum. We give an example in { IV 



II. 



SUPERHORIZON EVOLUTION WITH 
TRANSPORT EQUATIONS 



The "transport" approach is another implementation 
of the separate universe approximation [21 El |9] . We re- 
strict to canonical kinetic terms but allow an arbitrary 
number of fields, and set c = h = 1. The reduced Planck 
mass is Mp = (SttG)^^. The species of light scalar fields 
are indexed with Greek labels, and we make use of the 
slow-roll approximation throughout. 

In the standard picture, the fluctuations associated 
with increasing wavenumbers pass sequentially outside 
the horizon. After a few e-folds, the decaying mode of 
each fluctuation has died away and it is treated as a clas- 
sical object [TU] (although a detailed understanding of 
the quantum-to-classical transition is still lacking). To 
determine the subsequent evolution we invoke the sepa- 
rate universe approximation. Consider any large spatial 
region, smoothed on a scale larger than the horizon. The 
separate universe approximation asserts that, at any po- 
sition, the coarse-grained flelds evolve as they would in 
an unperturbed universe, up to gradient-suppressed cor- 
rections. 

The observable predictions of any inflationary model 
are not the fields themselves, but their correlation func- 
tions. In Refs. [HI in] these were studied using the distribu- 
tion of coarse-grained fields obtained after smoothing on 
the scale k^. To leading order, inflation predicts that this 
probability distribution is Gaussian and can therefore be 
characterized by the expectation value {4>a) = and 
variance {S(j>a^4'p) = Sq^ of the coarse-grained fields, 
where d(l)a = (pa ~ '^a- These one- and two-point func- 
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tions evolve according tc{^ 
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a 



(la) 
(lb) 



where we have uniformly omitted correlation functions 
of third-order or higher. The u-tensors are related to 
the field velocity (f)'^ and can be obtained using pertur- 
bation theory; for details, see Refs. O [S]. In general, 
each n-point function will be sourced by all other cor- 
relation functions, and for computational purposes the 
hierarchy must be truncated. In (la I, the truncation 



of three-point functions and above leaves a "loop cor- 
rection" to the one-point function, which is negligible in 
practical examples [101 US]. Eqs (la|-(lb| should be sup- 



plemented by boundary conditions obtained from match- 
ing to the known n-point functions computed perturba- 
tively at horizon-crossing. We will encounter an explicit 
example in §III| below. 

For comparison with experiment we require the cor- 
relation functions of the curvature perturbation, which 
seeds the observable density fluctuation [16]. Eqs. (lal~ 
( lb ) apply for any choice of correlation functions, but 
they are simplest to solve when our choice makes the sys- 
tem autonomous. A suitable choice is the n-point func- 
tions of the canonically normalized scalar fields whose 
potential energy supports the inflationary era. The cor- 
relation functions of the curvature perturbation, can 
be obtained from these by a gauge transformation. At 
fixed cosmic time t — tc, this transformation can be writ- 
ten 



^S(j)p - {S(j)aS(t>l3)) + ■ 



(2) 



where iV = J H dt measures the number of e-foldings of 
infiation and the omitted terms are higher-order in pow- 
ers of the small fluctuations All quantities in ([2| 
depend only on data local to the hypersurface t — t^.. 
Unless loop corrections are being retained, it is only nec- 
essary to keep contributions to the 2-point function from 



^ These equations were given for an inflationary model with an 
arbitrary number of fields in Refs. 0|9]. In the single-field case 
they appear in Ref. 1111 . 

After extinction of the decaying mode we are neglecting any 
further quantum effects. In principle, the probability distribution 
of the coarse-grained fields should be inherited from the gener- 
ating function of correlation functions, ie., the effective action. 
In a closed-time-path formulation there is no effective action in 
the usual sense, because its evolution (and that of the correlation 
functions) becomes stochastic. See, for example, Ref. 1121 . Over 
suffici ent ly lon g times this stochasticity cannot be neglected, and 
Eqs. ( |la[ |— ( jlbl l should be supplemented by Langevin terms, as 
originally proposed by Starobinsky 13 . In the full quantum 
theory, one would expect the evolution of correlation functions 
to be determined by a Schwinger-Dyson hierarchy. The transi- 
tion to a Boltzmann-Langevin hierarchy was studied by Calzetta 
& Hu [M]. 



N^a, dropping N^ap and any higher Taylor coefficients 
|15j . Also, N,a can be written in terms of the inflation- 
ary potential, V(4>), using 
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(3) 



Eq. (3| applies for an arbitrary V [6l IH [17] but invokes 
the slow-roll approximation. (The summation convention 
is suspended for ([s]), but used elsewhere in this paper.) 
The two-point function of C can be calculated using [2 [3] 



(CO - iV,a^,/3Sa/3. 



(4) 



Eqs. ([Ta|-([Tb]) and give a practical algorithm 

with which to study evolution of the power spectrum on 
the smoothing scale fc^.. By carrying out the computation 
for a few closely spaced kg one could numerically study 
the fcs dependence of the result. But this approach would 
be undesirable, requiring (la)-(lb) to be integrated sev- 
eral times for each model of interest. For example, if the 
resulting inflationary observables are to be used in draw- 
ing a sample from the space of models (perhaps as input 
to a code such as CosmoMC or MultiNest, which may 
require exponentially many samples), then the resulting 
time penalty can become appreciable. 



III. TRANSPORT OF THE SPECTRAL INDEX 



To do better, we return to Eq. ( lb I and use it to obtain 



a transport equation for the scale dependence directly. 
We write 

(C(ki)C(k2)) = (2^)H\V^ + k2)^ (5a) 



A„(ki)(50^(k2)) = (2^)^5^(ki+k2) 



2k\ 



(5b) 



Up to a normalization, the function P is the power spec- 
trum of the curvature perturbation, and can be computed 
using P = N ,aN ^pTjap- The fc-label argument of P and 
Sq^ should be associated with the smoothing scale ks 
discussed in Sjll] To study the variation of P with kg at 
equal times, it is conventional to use the parametrization 



P(fc,)~P,(fc,) 1 



l)*ln- 



where the spectral index Us satisfles 

dlnP(fc,.) 



1)* 



dln/c. 



(6) 



(7) 



k,=k. 



The "pivot" scale k^ is arbitrary and can be chosen for 
convenience. It is analogous to the arbitrary renormal- 
ization scale in scattering calculations. Observable quan- 
tities would be independent of our choice were all terms 
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in the expansion ([6| to be kept. In practice, after trun- 
cating at finite order it should be chosen approximately 
equal to the scale of interest, kg, to avoid an unwanted 
large logarithm. 

The derivative with respect to In kg is to be evaluated 
at equal times. Directly differentiating the 2-point func- 
tion of C, we obtain 



P din A:, 



(8) 



where we have introduced the matrix n^j^ = dS^^ /d In fc, 
which measures scale dependence of Eq,^. Note that 
the gauge-transformation factors N^a are A:s-independent. 
They depend only on the typical trajectory followed by 
the coarse-grained fields, and the time of evaluation. The 
same is true for Uap, meaning that it is simple to obtain 
an evolution equation for Uap, 



dn, 



d dE 



al3 



dN 



d\nk dN 



(9) 



This is of precisely the same form as the general transport 
equation (lb). Its solution will differ from E^^ through 



a different choice of boundary conditions. 

Eq. ([9]) is sufficient to determine the evolution of Ug 
from horizon crossing up to the end of inflation. It is 
only necessary to specify the initial value of ria^. To 
leading order in the slow-roll approximation, the fluc- 
tuations of canonically normalized scalar fields are un- 
coupled at horizon crossing, making E^^ proportional 
to a Kronecker-(5. With the normalization of (5b I, we 



have Eq,/3|* = H^Sap- Consider the pivot scale fc* and a 
nearby shorter mode with wavenumber k = fc* (1 + (5 In fc) 
and 6\nk > 0. When their decaying modes are lost, as 
described above, the fluctuations associated with these 
wavenumbers settle down to classical perturbations with 
amplitudes which differ by [TH] 



0.9 



P/H. 
60 

60 

40 

30 



FIG. 1. Superhorizon evolution of the spectral index (top) 
and power spectrum of curvature perturbations (bottom) for 
the double quadratic potential. 
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FIG. 2. Superhorizon evolution of the parameter /ml for the 
double quadratic potential. 



IV. EXAMPLES 
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S In fc. 



(10) 



These two modes cross the horizon at slightly different 
times, and therefore when compared at the same time 
the longer mode fc* experiences slightly more evolution. 
It follows that there is an extra displacement (SEq^ w 
— {dT,ap/dN)Slnk, because Sink measures the number 
of e-folds which elapse between horizon exit of the two 
modes. We conclude 



nal3\ 



(11) 



Note that Uap is symmetric when time is measured in 
e-folds, N. 

Eqs. are all that is required to compute the 

evolution of the spectral index. In the next section we 
present some illustrative examples, drawn for simplicity 
from the class of 2-field models. However, we emphasize 
that this method is valid for any number of scalar fields. 



Double quadratic inflation. Our first example is the 
well-studied model of double quadratic inflation [71 [TS] 



^22 



(12) 



where the mass ratio m^/m^ = 9 and m^^ = 10 ^Mp. 
We begin the calculation 50 e-folds before the end of in- 
flation, with (j) = 8.2Mp and x = 12.9Afp. In this model 
there is a single turn in field space. The turn is asso- 
ciated with a 'spike' in the amplitude of the local-mode 
bispectrum, typically measured by the parameter /nl; 
see Refs. [SI [H]. The turn is also associated with a 
step in the amplitude of the power spectrum, P, and a 
similar step in the spectral index. In Figs. [T] and [2| these 
behaviours can be compared. We leave a careful analysis 
for future work. 

Random potentials: nonmonotonic features. In 

this example we study an example of nonmonotonic evo- 
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FIG. 3. Superhorizon evolution of the spectral index (top) 
and power spectrum of curvature perturbations (bottom) for 
a particular realization of a randomly generated potential. 
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FIG. 4. Superhorizon evolution of the spectral index (top) 
and power spectrum of curvature perturbations (bottom) for 
the exponential quadratic potential. 



lution in the power spectrum. Nonmonotonic behaviour 
can occur for potentials which give rise to evolution less 
trivial than a single turn. We consider 

V ^ ^ [amn cos {mcj) + nx) + bmn sin (m0 + nx)\ ■ 

l^n^3 

(13) 

The coefficients amn and 5m„ are drawn randomly, mak- 
ing a low-order Taylor approximation to a "generic" po- 
tential. Frazer & Liddle [4] used potentials of this kind 
to study the generic features of inflation in a randomly 
generated landscape. 

Fig. [3] shows the evolution of the power spectrum and 
spectral index in a specific realization of this model, both 
of which are nonmonotonic. Unlike the simple case of 
double quadratic inflation, the spectral index does not 
mirror the evolution of the power spectrum. This is a 
consequence of the different initial conditions for riap 
and Sq^, which make n^p grow more slowly than Sq/j. 
Eq. ([s]) implies that the net result is a decrease in — 1. 

Exponential quadratic potential. Our final exam- 
ple exhibits a potential with unusual characteristics: the 
isocurvature modes do not decay and are present at the 
end of inflation. Therefore the curvature perturbation C, 
will continue to evolve: it has not yet reached its "adia- 
batic limit" [2D]. The observable statistical properties of 
C in models of this type will depend on unknown details 
of the post-inflationary history, including the details of 
reheating. Therefore such models are less predictive than 
those which converge to their adiabatic limit during in- 
flation. However, this is a question of calculability rather 
than principle. Models of this kind need not be intrinsi- 
cally less interesting from the point of view of fundamen- 
tal physics. 



We use a model introduced by Byrnes et al. [3T], and 
studied in more detail by EUiston et al. f5^. The poten- 
tial is 

V^V^X^e-^^\ (14) 

We choose Vq — 10~^Mp and A = 0.05, making our 
results comparable to Ref. [50] (the field values 55 e- 
folds before the end of inflation are (f) — O.OOlMp and 
X = 16Afp). Fig. [4] shows that the power grows strongly 
from roughly 45 e-folds after horizon crossing, reflecting 
the divergence of nearby trajectories in field space de- 
scribed by Elliston et al. The corresponding evolution 
of Us is a stretched step. The amplitude of P grows in 
comparison with NccNpriap, and the spectral index de- 
creases significantly. Its value at the end of inflation is 
observationally unacceptable, but nothing can be con- 
cluded until the post-inflationary history is specified. 

V. SUMMARY 

Inflationary models motivated by high-energy physics 
frequently involve many light scalar fields with compli- 
cated potentials. For this reason, it is important to find 
efficient numerical tools with which to predict observable 
parameters. 

In this paper, we have presented a simple method to 
compute the evolution of the spectral index in models 
where multi- field dynamics must be taken into account. 
Our approach uses a system of "transport equations" — a 
set of ordinary differential equations — to evolve the spec- 
tral indices nap of the 2-point functions of the canonically 
normalized scalar fields. These can be assembled into 
at the time of interest. 
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The transport equation for the indices Uap is exactly 
the same as the equation used to transport the power 
spectra of the field perturbations. However, owing to 
its different initial conditions, its evolution is not triv- 
ially related to the evolution of the power spectra. Its 
behaviour depends on the characteristics of the poten- 
tial. In pV| we exhibited specific examples that illustrate 
this behaviour. Although our examples were drawn from 
the class of 2-field models for reasons of simplicity, the 
method is valid for any number of fields. 

In a model with N fields, we require N integrations 
to track the evolution of the centroid ^c, and N{N + 
l)/2 to evolve the two-point function Eq^. To obtain 
a simple estimate, we assume the computational time 
is dominated by integration, and that integrating each 
component of Sq/j. Under these assumptions we can 
expect the time required to scale like N{N + 3)/2. 

Our method requires the introduction of another 
N{N -f l)/2 equations to track the evolution of Uap. As- 
suming integration of Uafj takes comparable time to Sa^g 
and <I>Q, (which should be a good approximation because 
and riai^ obey the same equation), the total time 
will now scale like N{N + 2). On the other hand, if we 
evolve the 2-point function at multiple scales in order to 
extract a numerical derivative by fitting a spline with m 



points, then we will require A'' -|- mN{N + l)/2 integra- 
tions. Typically m will be of order a few. 

These simple estimates show that, asymptotically, ei- 
ther method has complexity 0{N'^). To obtain reason- 
able accuracy when extracting a numerical derivative, 
one might wish to fit a spline with several points. In this 
case one might expect the transport approach to be more 
efficient by a modest numerical factor. Although not dra- 
matic, a small speed increase would certainly be worth- 
while when accumulated over a large number of samples. 

A similar approach could be adopted to calculate the 
scale dependence of any n-point function. The case of ob- 
servational interest is likely to be the bispectrum, where 
one might wish to compute the running of the local-mode 
amplitude /nl [22] • We leave this for future work. 
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